Replication

# dm and data visualization
library(tidyverse)

package 㤼㸱tidyverse㤼㸲 was built under R version 3.4.4-- Attaching packages --------------------------------------- tidyverse 1.2.1 --
v ggplot2 3.1.0     v purrr   0.2.5
v tibble  1.4.2     v dplyr   0.7.8
v tidyr   0.8.2     v stringr 1.3.1
v readr   1.2.1     v forcats 0.3.0
package 㤼㸱ggplot2㤼㸲 was built under R version 3.4.4package 㤼㸱tibble㤼㸲 was built under R version 3.4.4package 㤼㸱tidyr㤼㸲 was built under R version 3.4.4package 㤼㸱readr㤼㸲 was built under R version 3.4.4package 㤼㸱purrr㤼㸲 was built under R version 3.4.4package 㤼㸱dplyr㤼㸲 was built under R version 3.4.4package 㤼㸱stringr㤼㸲 was built under R version 3.4.4package 㤼㸱forcats㤼㸲 was built under R version 3.4.4-- Conflicts ------------------------------------------ tidyverse_conflicts() --
x purrr::accumulate() masks foreach::accumulate()
x tidyr::expand()     masks Matrix::expand()
x dplyr::filter()     masks stats::filter()
x dplyr::lag()        masks stats::lag()
x dplyr::select()     masks MASS::select()
x purrr::when()       masks foreach::when()

# dm and data visualization
library(tidyverse)
# bright color blind palette:  https://personal.sron.nl/~pault/ 
cblack <- "#000004FF"

data(washb_bang_tr)
data(washb_bangladesh_enrol)
data(washb_bangladesh_anthro)

washb_bang_tr = washb_bang_tr
washb_bangladesh_enrol = washb_bangladesh_enrol
washb_bangladesh_anthro = washb_bangladesh_anthro

latrines <- merged_bangladesh_data %>%
  select(block,laz,latown,latseal,covariates) %>%
  mutate(latrine = factor(ifelse(latown==0 | is.na(latown),"No latrine","Latrine")),
         improved_latrine = factor(ifelse(latown==1 & latseal==1 & !(is.na(latseal)),"Improved latrine","No improved latrine")))

latrines$latcat <- factor(NA, levels = c("No latrine","Latrine no water seal","Latrine with water seal"))
latrines$latcat[latrines$latrine=="No latrine"] <- "No latrine"
latrines$latcat[latrines$latrine=="Latrine" & latrines$improved_latrine=="No improved latrine"] <- "Latrine no water seal"
latrines$latcat[latrines$latrine=="Latrine" & latrines$improved_latrine=="Improved latrine"] <- "Latrine with water seal"
latrines

p1 <- ggplot(data=latrines, aes(x=laz,group=latcat,fill=latcat)) +
  geom_density(alpha=0.4) +
  scale_fill_manual(values=c(corange,cmagent, cgreen),
                    guide=guide_legend(title="Sanitation at enrollment")) +
  theme_minimal()
p1

# unadjusted
b_glmest_unadj_seal <- washb_glm(Y=latrines$laz,tr=latrines$latcat,id=latrines$block,
                            contrast=c("No latrine","Latrine with water seal"))

-----------------------------------------
 GLM Fit: Latrine with water seal vs. No latrine 
-----------------------------------------
                              Coef.      2.5%     97.5% Std. Error  z value          Pr(>|z|)
trLatrine with water seal 0.5200683 0.3440009 0.6961356  0.0898303 5.789453 0.000000007061599

# adjusted
b_glmest_adj_seal <- washb_glm(Y=latrines$laz,tr=latrines$latcat,W=latrines[covariates],id=latrines$block,
                          contrast=c("No latrine","Latrine with water seal"))

-----------------------------------------
Dropping 6 observations due to missing values in 1 or more variables
 Final sample size: 706 
-----------------------------------------

-----------------------------------------
Pre-screening the adjustment covariates:
-----------------------------------------

Likelihood Ratio Test P-values:
               P-value
month          0.83067
aged           0.76364
sex            0.79242
birthord       0.00000
momage         0.30098
momedu         0.00000
momheight      0.00000
Nlt18          0.06824
Ncomp          0.10751
watmin         0.64017
hfiacat        0.00013
elec           0.00000
floor          0.00000
walls          0.53965
roof           0.45929
asset_wardrobe 0.00000
asset_table    0.00103
asset_chair    0.00000
asset_khat     0.00000
asset_chouki   0.00243
asset_tv       0.00000
asset_refrig   0.00000
asset_bike     0.09432
asset_moto     0.00000
asset_sewmach  0.02966
asset_mobile   0.01257


Covariates selected (P<0.2):
                                       P-value
birthord       0.00000272283588528772020121750
momedu         0.00000060182286073082100316851
momheight      0.00000000000000000000001592703
Nlt18          0.06824241188032699467758135370
Ncomp          0.10751114688132597241310151048
hfiacat        0.00012873171546082565978169621
elec           0.00000016246202404920964655361
floor          0.00000000027912209753552275487
asset_wardrobe 0.00000000000042121291982584921
asset_table    0.00102929328603775654960017505
asset_chair    0.00000023697864269963536316445
asset_khat     0.00000000000012900393183580379
asset_chouki   0.00243206950221623368832424283
asset_tv       0.00000000003820491717153050757
asset_refrig   0.00000000081977785674527746812
asset_bike     0.09432420046593500972953449946
asset_moto     0.00000052972209002626561074809
asset_sewmach  0.02965887124644659908478949717
asset_mobile   0.01256587594867407051135721474


-----------------------------------------
 GLM Fit: Latrine with water seal vs. No latrine 
-----------------------------------------
                              Coef.        2.5%     97.5% Std. Error  z value   Pr(>|z|)
trLatrine with water seal 0.2001854 -0.00924158 0.4096123  0.1068505 1.873509 0.06099809

 Coef of covariates
                                      Coef.        2.5%       97.5%
birthord                        -0.09538543 -0.17217386 -0.01859699
momeduPrimary (1-5y)             0.00635140 -0.21197349  0.22467630
momeduSecondary (>5y)           -0.03971910 -0.24010021  0.16066201
momheight                        0.06447458  0.04986347  0.07908569
Nlt18                           -0.02794568 -0.09792990  0.04203853
Ncomp                           -0.00053605 -0.00995341  0.00888131
hfiacatMildly Food Insecure     -0.13362088 -0.40419275  0.13695100
hfiacatModerately Food Insecure  0.01066081 -0.17962502  0.20094664
hfiacatSeverely Food Insecure    0.07897930 -0.20969162  0.36765022
elec                             0.09195769 -0.07198243  0.25589782
floor                            0.19990013 -0.01302043  0.41282068
asset_wardrobe                   0.22448134  0.01420522  0.43475745
asset_table                     -0.00197888 -0.18968049  0.18572272
asset_chair                      0.10091079 -0.07456292  0.27638450
asset_khat                       0.20452258  0.03028987  0.37875528
asset_chouki                    -0.05429240 -0.24306253  0.13447772
asset_tv                         0.03289008 -0.15342942  0.21920959
asset_refrig                     0.16140511 -0.10846332  0.43127355
asset_bike                       0.02139492 -0.13438373  0.17717357
asset_moto                       0.09503119 -0.13120333  0.32126571
asset_sewmach                   -0.01352269 -0.24143107  0.21438569
asset_mobile                    -0.03592414 -0.23638673  0.16453846

# unadjusted
b_glmest_unadj_no_seal <- washb_glm(Y=latrines$laz,tr=latrines$latcat,id=latrines$block,
                            contrast=c("No latrine","Latrine no water seal"))

-----------------------------------------
 GLM Fit: Latrine no water seal vs. No latrine 
-----------------------------------------
                            Coef.        2.5%     97.5% Std. Error   z value  Pr(>|z|)
trLatrine no water seal 0.0382343 -0.08399555 0.1604641 0.06236217 0.6131008 0.5398097

b_glmest_adj_no_seal <- washb_glm(Y=latrines$laz,tr=latrines$latcat,W=latrines[covariates],id=latrines$block,
                          contrast=c("No latrine","Latrine no water seal"))

-----------------------------------------
Dropping 4 observations due to missing values in 1 or more variables
 Final sample size: 900 
-----------------------------------------

-----------------------------------------
Pre-screening the adjustment covariates:
-----------------------------------------

Likelihood Ratio Test P-values:
               P-value
month          0.21574
aged           0.86874
sex            0.85603
birthord       0.00007
momage         0.05102
momedu         0.00000
momheight      0.00000
Nlt18          0.02342
Ncomp          0.78382
watmin         0.58348
hfiacat        0.01200
elec           0.00021
floor          0.00002
walls          0.67221
roof           0.54869
asset_wardrobe 0.00000
asset_table    0.18722
asset_chair    0.00243
asset_khat     0.00000
asset_chouki   0.00122
asset_tv       0.00000
asset_refrig   0.00061
asset_bike     0.03828
asset_moto     0.03734
asset_sewmach  0.09164
asset_mobile   0.02607


Covariates selected (P<0.2):
                    P-value
birthord       7.129668e-05
momage         5.101940e-02
momedu         7.562345e-07
momheight      9.572032e-32
Nlt18          2.342201e-02
hfiacat        1.200262e-02
elec           2.120838e-04
floor          1.579166e-05
asset_wardrobe 1.193311e-06
asset_table    1.872224e-01
asset_chair    2.431295e-03
asset_khat     1.793530e-06
asset_chouki   1.224883e-03
asset_tv       9.867935e-07
asset_refrig   6.068402e-04
asset_bike     3.827750e-02
asset_moto     3.734111e-02
asset_sewmach  9.164459e-02
asset_mobile   2.606618e-02


-----------------------------------------
 GLM Fit: Latrine no water seal vs. No latrine 
-----------------------------------------
                              Coef.       2.5%      97.5% Std. Error    z value  Pr(>|z|)
trLatrine no water seal -0.03123738 -0.1532639 0.09078917 0.06225844 -0.5017372 0.6158524

 Coef of covariates
                                      Coef.        2.5%      97.5%
birthord                        -0.06458955 -0.13292813 0.00374903
momage                           0.00447556 -0.01227145 0.02122257
momeduPrimary (1-5y)            -0.02593891 -0.21883467 0.16695685
momeduSecondary (>5y)            0.01549283 -0.16744499 0.19843066
momheight                        0.06601027  0.05511840 0.07690213
Nlt18                           -0.01208007 -0.06943355 0.04527341
hfiacatMildly Food Insecure     -0.01719910 -0.22748631 0.19308811
hfiacatModerately Food Insecure -0.02618838 -0.18497808 0.13260133
hfiacatSeverely Food Insecure   -0.10771005 -0.39303945 0.17761934
elec                             0.01610124 -0.11062790 0.14283039
floor                            0.27939377  0.08078173 0.47800581
asset_wardrobe                   0.20580959  0.01544444 0.39617474
asset_table                     -0.07516113 -0.23042631 0.08010405
asset_chair                      0.06931524 -0.09129706 0.22992754
asset_khat                       0.05792481 -0.10367188 0.21952150
asset_chouki                    -0.07299354 -0.23527860 0.08929152
asset_tv                         0.11700790 -0.03708876 0.27110456
asset_refrig                     0.13955711 -0.16989073 0.44900496
asset_bike                       0.10574671 -0.05245537 0.26394879
asset_moto                      -0.01579300 -0.25348966 0.22190367
asset_sewmach                    0.11350529 -0.16726049 0.39427106
asset_mobile                     0.03609640 -0.13517109 0.20736389

unlist(b_tmleest_adj_seal$estimates$ATE)

        psi     var.psi         CI1         CI2      pvalue 
0.213195358 0.010763582 0.009849898 0.416540818 0.039884714 

unlist(b_tmleest_adj_no_seal$estimates$ATE)

        psi     var.psi         CI1         CI2      pvalue 
-0.02888098  0.00494765 -0.16674647  0.10898451  0.68137042 

Kenya

nrow(merged_kenya_data)

[1] 8928

# unadjusted
k_glmest_unadj <- washb_glm(Y = kenya_data$laz,tr=kenya_data$latcat,id=kenya_data$block,
                            contrast=c("No improved latrine","Improved latrine"))

-----------------------------------------
 GLM Fit: Improved latrine vs. No improved latrine 
-----------------------------------------
                       Coef.      2.5%     97.5% Std. Error  z value      Pr(>|z|)
trImproved latrine 0.2479052 0.1232202 0.3725903 0.06361481 3.896974 0.00009740218

# adjusted
k_glmest_adj <- washb_glm(Y=kenya_data$laz,tr=kenya_data$latcat,W=kenya_data[kenya_covariates],id=kenya_data$block,
                          contrast=c("No improved latrine","Improved latrine"))

-----------------------------------------
Dropping 164 observations due to missing values in 1 or more variables
 Final sample size: 1937 
-----------------------------------------

-----------------------------------------
Pre-screening the adjustment covariates:
-----------------------------------------

Likelihood Ratio Test P-values:
            P-value
month       0.41581
HHS         0.01380
agey        0.00602
sex         0.00001
mother_age  0.00017
motherht    0.00000
mother_edu  0.00000
u18         0.00239
Ncomp       0.04760
water_time  0.59419
floor       0.00724
roof        0.14014
cow         0.00625
goat        0.13611
chicken     0.00988
dog         0.18957
elec        0.00649
radio       0.16068
tv          0.00103
mobilephone 0.00135
clock       0.02504
bicycle     0.31627
motorcycle  0.58492
stove       0.00000


Covariates selected (P<0.2):
                 P-value
HHS         1.379557e-02
agey        6.016331e-03
sex         1.490603e-05
mother_age  1.714750e-04
motherht    3.887960e-31
mother_edu  8.415038e-09
u18         2.387220e-03
Ncomp       4.760270e-02
floor       7.244095e-03
roof        1.401401e-01
cow         6.249603e-03
goat        1.361061e-01
chicken     9.881421e-03
dog         1.895659e-01
elec        6.493348e-03
radio       1.606750e-01
tv          1.029531e-03
mobilephone 1.349252e-03
clock       2.503707e-02
stove       1.520519e-06


-----------------------------------------
 GLM Fit: Improved latrine vs. No improved latrine 
-----------------------------------------
                       Coef.       2.5%     97.5% Std. Error  z value    Pr(>|z|)
trImproved latrine 0.1619255 0.04552707 0.2783239 0.05938695 2.726617 0.006398718

 Coef of covariates
                                   Coef.        2.5%       97.5%
HHSModerate                  -0.17349766 -0.32784410 -0.01915123
HHSSevere                     0.26122019 -0.16152382  0.68396421
agey                         -0.42676567 -0.66907277 -0.18445857
sexMale                      -0.21970055 -0.31158411 -0.12781698
mother_age                    0.02197288  0.01285362  0.03109214
motherht                      0.04004494  0.02831249  0.05177739
mother_eduAny secondary       0.18679715  0.07503326  0.29856104
mother_eduComplete primary    0.06635516 -0.04110909  0.17381941
u18                          -0.06228162 -0.09267246 -0.03189078
Ncomp                        -0.00451917 -0.01299730  0.00395895
floorCement/concrete         -0.00694482 -0.21659605  0.20270640
floorMissing/DK              -1.00035315 -1.27069155 -0.73001475
roofIron/other               -0.00036013 -0.10892189  0.10820162
cow                           0.02128260 -0.01563309  0.05819829
goat                          0.01228318 -0.02937551  0.05394187
chicken                       0.00183331 -0.00334463  0.00701126
dog                           0.01876128 -0.03584364  0.07336620
elecHas electricity           0.09812476 -0.10920199  0.30545151
radioHas radio               -0.00708682 -0.11170660  0.09753296
tvOwns TV                     0.10017201 -0.02940975  0.22975376
mobilephoneAny mobile phones  0.10436950 -0.04396080  0.25269979
clockHas clock               -0.00749205 -0.14445876  0.12947466
stoveHas stove                0.12874814  0.01422996  0.24326632

unlist(k_tmleest_adj$estimates$ATE)

        psi     var.psi         CI1         CI2      pvalue 
0.145983581 0.004853697 0.009433355 0.282533807 0.036135141 

tbl

                           Population (%) Mean LAZ (SD)    Unadjusted Difference in LAZ (95% CI) P-value   Adjusted Difference in LAZ - GLM (95%) P-value 
No latrine                 "513 (47%)"    "-1.8939 (0.98)" "Ref"                                 "Ref"     "Ref"                                  "Ref"   
Latrine with no water seal "391 (35%)"    "1.8557 (1.00)"  "0.0382 (-0.08, 0.16)"                "0.5398"  "-0.0312 (-0.15, 0.09)"                "0.6154"
Latrine with water seal    "199 (18%)"    "-1.3739 (1.01)" "-0.5200 (0.34, 0.69)"                "<0.0001" "0.2002 (-0.06, 0.41)"                 "0.610" 
No improved latrine        "1737 (83%)"   "-1.578 (1.08)"  "Ref"                                 "Ref"     "Ref"                                  "Ref"   
Improved latrine           "364 (17%)"    "1.3301 (1.08)"  "0.2479 (0.12, 0.37)"                 "<0.0001" "0.1619 (0.05, 0.28)"                  "0.0064"
                           Adjusted Difference in LAZ - TMLE (95%) P-value 
No latrine                 "Ref"                                   "Ref"   
Latrine with no water seal "-0.0289 (-0.17, 0.12)"                 "0.6814"
Latrine with water seal    "-0.2132 (0.01, 0.42)"                  "0.0399"
No improved latrine        "Ref"                                   "Ref"   
Improved latrine           "0.1461 (0.01, 0.29)"                   "0.0361"